Spatial fluctuations of the pairing potential in 
disordered superconductors 



Grzegorz Litak 

Department of Mechanics, Technical University of Lublin Nadbystrzycka 36, 

PL-20-618 Lublin, Poland. 



Abstract 

We study the effect of site diagonal, non-magnetic, disorder on tlie a pairing ampli- 
tude in an extended Hubbard model with the intersite attraction. Analyzing fluc- 
tuations of a pairing potential we discuss the instability of mixed solutions, 's + d' 
and 's + \d\ in presence of disorder. The influence of disorder on extended s- and 
d-wave superconductors appear to be comparable but in certain regions of the phase 
diagram, even weak disorder can change the symmetry of the order parameter. 
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1 Introduction 

The study of disordered superconductors witli exotic pairing is of general 
interest. Z)-wave pairing plays an important role in the superconductivity of 
the cuprates [1] while p-wave in strontium ruthenate (Sr2Ru04) [2]. Depending 
on symmetry of order parameter, the effect of disorder is dramatically different 
on the superconducting states which are described as having isotropic, s-wave, 
or anisotropic, extended s-, d- and p-wave, order parameter symmetry [3-13]. 
The works on disorder effect on superconductors are stimulated by the recent 
experiments [14-19]. 

In this paper we analyze an extended Hubbard model with intersite attrac- 
tion, whose phase diagram includes both s- and d-wave regions and introduce 
disorder into the problem by allowing the site energies, £j, to be independent 
random variables. 
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In the case of anisotropic s- and d-wave pairing there is no strict equivalent of 
the Anderson Theorem [3,5,9] which governs the response of isotropic s-wave 
superconductors to randomness in the crystal potential. However, a weaker 
statement, which may be similarly useful, can be formulated [7,11,12]. Here 
we will calculate the ensemble average < A^A*; > where Ajj is the pairing 
amplitude in the cases where the sites, i,j and I, are nearest neighbours to 
lowest order in the fluctuations of the site energies Ssi = Ei — Eq about their 
mean Eq. Our results imply that the disorder does not induce large fluctuations 
in the pairing potentials and hence for systems of large coherence length ^ the 
amplitude of pairing potential |A(ij)| may rather to be the same for all bonds 
in X and y directions [9,11,12]. Eventually, we compare the corresponding 
standard square deviation of fluctuating potentials cr{6Aij) and (Fig. 
1). Their ratio T = {a{6Aij)/a{5Ei)y will be a criterion of pairing potential 
fluctuations leading flnally to the destruction of superconducting phase. We 
shall be also interested what is the effect of disorder on anisotropic s- and d- 
wave superconductors, where they coexist, and if the disorder favour any of 
particular solution. 




6e , 6A 



Fig. 1. The examples of Gaussian distribution for Ssi with a mean square deviation 
a{Sei) = < > and the expected distribution of a pairing potential (5Aj with 

corresponding a{5Aij) = ^< |Ajjp >. Fluctuating potentials Ssi and 5Ajj are in 

the units of \ < e? > 



The paper is organized as follows. In Sec. 2 we present the model in the clean 
limit, introduce the Hamiltonian and approximations used in paper. Here we 
also discuss briefly the phase diagram of the clean system. In Sec. 3 we include 
disorder induced spatial fluctuations of a pairing potential. Here we discuss 
the formalism used in our analysis and show the results on the disordered 
effect on anisotropic superconductors. Finally we investigate the stability of 
various solutions in presence of disorder. Section 4 contains conclusions and 
remarks. 
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2 The model in clean limit. 



Here we use a single band, extended, Hubbard model with an intersite attrac- 
tion. It is defined by the Hamiltonian [20] : 

H = + ^ijHaCja + Yl ^^ia^ja' ' I^Yl 4aCia, (1) 

ijcT ijacr' i<7 



where and are usual, fermionic operators which create and annihilate, 
respectively, an electron with spin a at the lattice site labeled by i, tij is a 
electron hopping integral, rii^ is the operator of particle number of spin a at 
site i, £i is the site energy, varying from site to site in random fashion, at the 
site i with mean value Sq —< Si >— Q and Wij is the interaction potential of 
two electrons with opposite spins on neighbour sites Finally, denotes 
the chemical potential. 

Hartree-Fock-Gorkov equation of motion for Gorkov Greens 2x2 functions 
G{i,j\uj) yields: 



E 



{uj - ei + fi)6ii - til - Wii{< cf^cii > -m) Wii < Ciici^ > 
Wii < Ci^cii > {uu + ei- + til + Wii{< Cj^Qi > -ni) 
xG{l,j;uj) = l6ij, (2) 



where we assumed a paramagnetic state riii = 



The charge rii on the site i and pairing potential An — Wu < Qj^c^ > for 
neighbour sites i, I can be expressed by following equations: 



rii 



Ail 



2 f 1 
/ do; Im Giifi, i: cu)— 5 , 

— OO 

W' ^ 1 

/ d^I^'^l2(^,/;a;)^^;^, (3) 



where /3 = and ks is the Boltzmann constant. 



3 



For a pure system we take the lattice Fourier transform of Eq. 2 with Si — 
and find: 



(4) 



where /x' is shifted chemical potential /j,' = /j, + Wn, and the Hartree-Fock 
kinetic energy efc is given by 

^k^ek- —2^nk-fk, (5) 
k 



where 



rik =< CfeC+ >= E < 4cj > e^(«^-^)'=. (6) 

ij 



Following Ref. [20] we will neglect the Fock term ^X]fc^fc7fc (Eq- 5). Thus, 
for two dimensional lattice and the electron hopping defined between nearest 
neighbour sites only, we get 

4«^;. = E^^ie'^''^"'^^' = -^7fe, (7) 



where 

7^, — 2(cosA;j; + cosky). (8) 

As can be readily shown the order parameter satisfies the following gap 
equation: 

A. = ^E^^tanh(^-j, (9) 

where 

iyfe_,= |iy|7fe_, (10) 
= \W\ {l^k^nj^ + 2sin(fc,)sin(g,) + 2sin(fc,)sin(g,)) 

and 

77fc = 2{coskx — cosky), (11) 
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while Eq denotes quasi-particle energy: 



Sq = ^eJ-A^ , ~eq = eq-i^'. (12) 

Generally, the singlet type of solution (Eq. 9) can be written as 

Ak^A^^^ + A'rj^ (13) 

for the real type solution if s and d parts of A (A;) have the same phase and 

Afe = A^7fe + iA'^77fe (14) 

for the complex solution. Pairing amplitudes: A*, A'' are defined as real num- 
bers, corresponding to s- and d-wave components respectively. 

The pairing parameters A^, A*^ (Eqs. 13, 14) may be determined from Eq. (9). 
Namely: 

where and a = 1 for a real solution {s + d) while a = i for a complex one 
{s + id). 

The above set of equations for pairing potentials (15) should be completed by 
the equation for the chemical potential /i (Eq. 3). It can be written as: 

n-l = -|-V^ tanhf^V (16) 



The gap equations (Eqs. 15) can be used to calculate the superconducting 
critical temperature Tc after the linearization. The two separate equations 
read: 



N^Si, [ 2 )' 
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where Pc^j^- 



In the Fig. 2 we present Tc resuhs for two dimensional lattice versus band 
filling n. The interaction parameter W was chosen as W/D — —0.3, —0.5, 
—0.7, where D — 8t denotes a bandwidth. Because of the assumed form of 




Fig. 2. The superconducting critical temperature Tc versus band filling n for s- 
and d-wave pairing (dashed and full lines respectively) for three values of intersite 
attraction \W\/D = 0.3, 0.5, 0.7 

a dispersion relation 6^ (Eq. 7) the model (Eq. 1) possesses the particle- 
hole symmetry, therefore Tc for band-filling n and 2 — n is the same. The 
results show that s-wave type of superconductivity exists for low electron or 
hole concentration while d-wave type exists close to half-band filling [20]. 
For a large enough electron interaction parameter W/D > 0.2 the curves 
Tc — Tc{n) for s- and d-wave pairing cross and regions of s- and d- wave 
superconductivity arc not separated. The pairing parameters, and A*^ at 
zero temperature (T = K) should be calculated from the full set of equations 
Eqs. (15). For the real type solution (Eq. 13) they are plotted for the same 
three values of the electron attraction W/D — —0.3, —0.5, —0.7 in Fig. 3a. 
Note that in all cases there are small but visible regions of mixed 's + d' 
solutions. In Fig. 3b we have plotted the same for the complex solution (Eq. 
14) of mixed 's + id' pairing. Clearly, the regions with possible mixed solution 
in case of 's + id' pairing are much larger than for 's + d' mixture. We will 
consider these two solutions more carefully. 

In the case of mixed solution, the order parameter (complex or real) in the 
lattice real space Ajj is defined as the sum of Af^- and A^^-. Namely: 

A,, = A^^. + a(-l)'-^Af„ (18) 



where / denotes the neighbour site (/ = j, j + 1, j + 2, j + 3 Fig. 4). a is equal 
to 1 and i for a real and complex solution, respectively. 
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Fig. 3. Amplitude of the pairing potential A versus band filling n for s- and d- wave 
pairing (dashed and full lines respectively) for three values of intersite attraction 
|VF|/Z) = 0.3, 0.5, 0.7. Note the difference around n = 0.4, figure (a) corresponds 
to the real solution 's + d' while (b) to the complex one 's + id' respectively. 
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Fig. 4. Neighbour sites on two a dimensional square lattice i, j and I = j + 1, j + 1, 
j + 2, J + 3. 

Prom the above equation (Eq. 18) one can determine (Fig. 4) the s and d 
components: 
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-(A,,- + A,,-+0 



(19) 



The free energy F for a finite temperature T can be calculated from the 
following formula [20]: 



[n 



l)ek - IksT \n [2 cosh 



2kuT 



|Afe(T)p 
W 



(20) 
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Fig. 5. Free energies F of mixed solutions: (a) the real combination 's + d' and 
(b) the complex one 's + id' for the same band filling n = 0.4 and the attraction 
parameter |W|/-D = 0.7. Additional arrows show the minima of free energies F. 

the corresponding derivatives determine the integral equations (Eqs. 14 and 
15): 

dF dF 



In Fig. 4 wc have plotted the free energies F of both mixed solutions: the real 
combination 's + d' (Fig. 5a) as well as the complex one 's + ici' (Fig. 5b) 
versus s- and d-wave amphtudes A* and A'*. The band filling n was chosen 
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to be n = 0.4 and the attraction parameter = 0.7 is large enough to 

produce the mixed sohitions (Figs. 2-3). Additional arrows in the plots show 
the minima of free energies F. Clearly they show the mixed type of solutions. 
To check which solution is more favorable we show also the difference between 
free energies for two superconducting mixed states AF = F(^s-\-d') — F(^s-\-\d') 
versus s- and d-wave amplitudes (Fig. 5c). One can see that in the region of 
coexistence AF is positive. Thus, the solution with s-\-\d symmetry is favored 
by smaller free energy F. 



3 Disorder induced Fluctuations of pairing potential 

To go further we apply the same strategy as in Refs. [9,10,21,22] we treat 
random site energies as perturbations and we proceed by solving the Dyson 
equation for the Gorkov matrix of the pure superconductor evaluated at a 
frequency uj 

G{t,r,uj) = G%t,r,u;) + J2G'{t,r,u)ViG{l,j;oo), (22) 

I 

where V; is the impurity potential matrix: 

-(::) 

To the lowest order in £j we get: 

G(z, j; uj) = G\i, j; ^) + ^ G°(z, n; uj)VnG\n, j; u). (24) 

n 

Following Eqs. (3) we express the interesting quantities Aj^ and An, where 
I and j are nearest neighbours of i, in the lowest order of Si perturbations 
by means of the disordered Green function (Eq. 24) and calculate the mean 
square deviation of the pairing parameter as: 

< SAijSAl >=< AijAl >-<Aij>< A* >, (25) 

where bonds ij and il in the pairing potentials A^j and An, denote bonds which 
can be chosen as parallel (then j and I coincide) as well as perpendicular. 

The assumption is that random site energies Si in Eqs. (1,22,23,24) are inde- 
pendent variables; then averages < > ,< £j > being independent of the site 
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index i and < EiEj >—< ef > 5ij lead to 



< 6\,6A;i >= < e] >, (26) 



where / = j for parallel {r{j — vjj) and I — j + 1 perpendicular bonds 
— r^), respectively (Fig. 4). Thus, fluctuations are determined by the 

corresponding coefficients rjj and F^. They can be calculated from Eqs. (3) 
and (22-26): 



X e 



X e 



w _ WuW,, 1 ^ 7 , Im(G°(fc;c.)T3G°(fc-g;c.))^^ 
y \ fc -oo 

jk' (Ri-Rj)^ ^iq(Rj-Ri) ^ (27) 
After evaluation at T = K it yields: 



'-^J AT ^ 



J du;lm{GUk;u;)GUk-q;u;) 
k -oo 

GUk;u;)Gl,{k-q-u;)}e-'''i^^-^^ 



5: /do;' Im {G',,(k';u;')GUk'-q;u;') 
k' -oo 

- G\^{k';u)Gl^{k' - q;u)}^^^^'-^^)\ M^^-^^). 



(28) 



Ntt 



Depending on site / we get the formula for Tjj if sites j and / arc identical {Ri 
=Rj) and F^- for / = j + 1 (Fig. 4). After the integration over u, and lu' (Eq. 
28) we get following formulae: 



F^. 



TV ^ 



2N 



y-v Afcefc + Afcefc _ q ^ 



^ik{Ri-R^) 



X 



2N^ 



-j^ {Ek + Ek- q)EkEk - q 



ik{R,-Ri) 



}Q(R3-Ri 



(29) 
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Fig. 6. Fluctuation parameter F = F^, versus band filling n for a pairing potential 
Aij for three values of intersite attraction |VF|/Z) = 0.3, 0.5, 0.7. Full lines corre- 



spond to the real solution 's + d' 
respectively. 



while dashed lines to the complex one 's + id' 



1 

N 



E 

q 



2N 



y Afcefc + Afcefc-q ^k(Ri-Ri) 

{Ek + -E'fe - q)EkEk - q 



In Fig. 6 we plot the coefficient versus band filling n. Full and dashed lines 
correspond to 's + rf' and 's + irf' solutions, respectively. One can see that in 
case of weak interactions the results are very similar, while for large interaction 
{W/D — —0.7) where mixed solutions are present (n ^ 0.4) we observe a 
large difference. Clearly, in case of a purely real solution the fluctuations are 
larger. Interestingly, we find PH = for half filled band n = 1. This result 
corresponds to the similar one in case of attractive 'on site' (negative U) 
interaction [9,10,21,22] and this is due to the particle hole symmetry. 

It is also possible to analyze fluctuations of each of the components of the 
superconducting order parameter (Eq. 17). The corresponding coefficients F^^- 
and rfj determine fluctuations of the pairing parameter of s-wave symmetry 
A|- and d-wave one Af-. Namely: 



< 



5A^5A^ 



F' 



< ej > 



< et > 



(30) 



where we have defined 



o(4 + ri^) 



(31) 
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Fig. 7. Fluctuation parameter T = F*, F versus band filling n for s- and d- 
wave pairing potential (full and dashed lines respectively) (dashed and full lines 
respectively) for three values of intersite attraction |VF|/D = 0.3, 0.5, 0.7. Figure 
(a) corresponds to the real solution 's + d' while (b) to the complex one 's + id' 
respectively. 

The parameters r*^'^) governing disorder induced fluctuations of the order pa- 
rameter ^s{d) arc plotted versus band filling in Fig. 7. Figure 7a corresponds 
to the real combination of 's + d' solution while Fig. 6b corresponds to the 
complex one 's + id\ Note that generally, fluctuations of 's + d' and d -\- id 
components are relatively small |F| < 0.012: F"* and F'' are of the same or- 
der. However their relative value IF^I/IF*^! depends on band flUing n, showing 
that the change of symmetry is possible by disorder. In case of particle-hole 
symmetry at half fllling both F* and F*^ go to zero (like F^ in Fig. 6). 



4 Conclusions and Discussion 

Solving the extended Hubbard model (1), for the appropriate system param- 
eters, we have analyzed possible singlet solutions with s- and rf-wavc order 
parameter. Especially we have concentrated on the existence of mixed, 's + d' 
and ^s-\-\d\ solutions. On account of that the Anderson theorem could not be 
applied to anisotropic superconductors we analyzed the disorder induced, spa- 
tial fluctuations of order parameter amplitudes A'' and A'^. Our results show 
that the disorder does not induce large fluctuations in the pairing potentials 
and the amplitude of pairing potential |A(ij)| may rather be the same for all 
bonds in x and y directions. However, in regions of mixed solution the fluctua- 
tions are larger for the solution, which imphes that 's-l-id' is more stable 
in the presence of disorder. Interestingly, we observed zero fluctuations limit in 
case of particle-hole symmetry as in the case of the negative U Hubbard model 
[9,10,21,22]. This result is obtained here at zero temperature (T = K) in the 
limit of very weak disorder. Nevertheless, for stronger non-magnetic disorder 
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we observe a pair breaking effect strongly influencing the critical temperature 
[11,12]. Moreover the effect is stronger if the chemical potential passes a Van 
Hove singularity. In our case it is a half-filled situation n = 1. 

Recently Ghosal et al. [23] have done self-consistent calculations of the pairing 
amplitude for finite cluster size and various strength of disorder. Their results, 
obtained away from half filling n ^ 1, show that in the case of strong enough 
disorder, large fluctuations of Ajj can lead to superconducting islands. In 
that situation the disorder can lead to phase fluctuations and to pseudogap 
phenomena [24] . To investigate these effects we need to use more sophisticated 
methods i.e. to go beyond the Hartree-Fock-Gorkov approximation as in Ref. 
[25]. 
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